Create a table summarizing data layers within each MTBS fire. Layers
were assembled and reduced in Google Earth Engine, then the summary file
exported to Google Drive. Downloaded the summary table locally from
Google Drive, then ran this script to visualize and explore the
data.
# how many fires in this data set?
message("# of fires in dataset: ")
## # of fires in dataset:
message(nrow(fire_df))
## 37
# create output directory for figures
if(!dir.exists(here::here("figures"))){
dir.create(here::here("figures"))
}
# create a table to summarize attributes of interest per fire
fire_table <- as.data.frame(fire_df) %>%
# select columns of interest
dplyr::select(Fire_Name, Year, Acres, gedi_coverage, lodgepole, ponderosa, spruceFir, disturbed_burned, disturbed_unspecific,
#disturbed_logged, regenerating_disturbed, regenerating_harvested
) %>%
# reorder the rows based on multiple variables
dplyr::arrange(desc(gedi_coverage), desc(lodgepole), desc(ponderosa), desc(spruceFir), desc(disturbed_burned)) %>%
# rename the columns for interpretation
dplyr::rename("Fire Name" = Fire_Name,
"Year" = Year,
"Acres" = Acres,
"GEDI %" = gedi_coverage,
"Lodgepole %" = lodgepole,
"Ponderosa %" = ponderosa,
"Spruce Fir %" = spruceFir,
"Disturbed Burned %" = disturbed_burned,
"Disturbed Unspecified %" = disturbed_unspecific,
#"Disturbed Logged %" = disturbed_logged,
#"Regenerating Disturbed %" = regenerating_disturbed,
#"Regenerating Harvested %" = regenerating_harvested
)
# add a categorical column with dominant forest type
temp <- fire_df %>%
dplyr::select(Fire_Name, lodgepole, ponderosa, spruceFir) %>%
tidyr::gather(major_forest_type, forest_percent, lodgepole:spruceFir) %>%
dplyr::group_by(Fire_Name) %>%
dplyr::slice(which.max(forest_percent)) %>%
dplyr::select(Fire_Name, major_forest_type) %>%
st_set_geometry(NULL)
fire_df <- fire_df %>%
dplyr::left_join(temp)
remove(temp)
kableExtra::kable(fire_table) %>%
kableExtra::kable_styling(bootstrap_options = "striped",
full_width = F,
position = "left")
|
Fire Name
|
Year
|
Acres
|
GEDI %
|
Lodgepole %
|
Ponderosa %
|
Spruce Fir %
|
Disturbed Burned %
|
Disturbed Unspecified %
|
|
MAD CREEK
|
2001
|
1435
|
99.7
|
1.1
|
0.0
|
74.6
|
0
|
0
|
|
BIG ELK
|
2002
|
4347
|
98.5
|
19.3
|
47.9
|
0.0
|
0
|
0
|
|
CRYSTALFIRE
|
2011
|
2855
|
95.3
|
2.7
|
54.2
|
0.0
|
0
|
0
|
|
LOWER NORTH FORK FIRE
|
2012
|
3436
|
93.8
|
0.2
|
31.2
|
0.0
|
0
|
0
|
|
#6
|
1989
|
1547
|
93.2
|
0.0
|
6.7
|
0.0
|
0
|
0
|
|
SNAKING
|
2002
|
2080
|
91.4
|
0.2
|
65.4
|
2.7
|
0
|
0
|
|
FOURMILE CANYON
|
2010
|
5865
|
91.1
|
1.2
|
67.0
|
0.0
|
0
|
0
|
|
BIG RED
|
2017
|
2868
|
90.0
|
80.6
|
0.0
|
0.9
|
0
|
0
|
|
UTE CREEK
|
1994
|
3225
|
89.2
|
2.1
|
0.0
|
36.0
|
0
|
0
|
|
PICNIC ROCK
|
2004
|
9014
|
88.2
|
0.0
|
15.9
|
0.0
|
0
|
0
|
|
MT. ZIRKEL COMPLEX (HINMAN)
|
2002
|
17627
|
86.4
|
10.1
|
0.1
|
23.1
|
0
|
0
|
|
BUFFALO CREEK
|
1996
|
9792
|
85.5
|
0.0
|
9.4
|
0.0
|
0
|
0
|
|
OLD MAN
|
2001
|
1169
|
83.4
|
0.0
|
0.0
|
0.0
|
0
|
0
|
|
MEADOW CREEK
|
2010
|
1432
|
83.3
|
3.8
|
0.1
|
8.9
|
0
|
0
|
|
MT. ZIRKEL COMPLEX (BURN RIDGE)
|
2002
|
15478
|
71.8
|
4.4
|
0.0
|
58.0
|
0
|
0
|
|
BOBCAT
|
2000
|
9131
|
69.8
|
1.4
|
75.9
|
0.0
|
0
|
0
|
|
BIG FISH
|
2002
|
17273
|
68.3
|
3.2
|
0.0
|
58.2
|
0
|
0
|
|
HAYDEN PASS
|
2016
|
17978
|
61.7
|
2.0
|
21.5
|
15.0
|
0
|
0
|
|
BEAVER CREEK
|
2016
|
44222
|
60.4
|
47.7
|
0.0
|
3.6
|
0
|
0
|
|
WALDO CANYON
|
2012
|
20112
|
58.6
|
0.3
|
47.2
|
0.2
|
0
|
0
|
|
LOST SOLAR
|
2016
|
5142
|
56.3
|
3.1
|
0.0
|
52.9
|
0
|
0
|
|
SUNNYSIDE
|
1989
|
1416
|
47.8
|
0.0
|
27.5
|
0.0
|
0
|
0
|
|
CANYON
|
1988
|
1092
|
47.5
|
0.0
|
40.8
|
0.0
|
0
|
0
|
|
HIGH PARK
|
2012
|
90769
|
46.8
|
13.9
|
49.6
|
1.1
|
0
|
0
|
|
LOST LAKES
|
2002
|
5888
|
44.9
|
0.2
|
0.0
|
51.8
|
0
|
0
|
|
UNNAMED
|
1989
|
1305
|
43.6
|
0.0
|
43.7
|
0.0
|
0
|
0
|
|
DUCKETT
|
2011
|
4705
|
41.2
|
0.4
|
17.4
|
0.0
|
0
|
0
|
|
COAL SEAM
|
2002
|
12004
|
37.3
|
0.1
|
0.0
|
0.5
|
0
|
0
|
|
OVERLAND
|
2003
|
3232
|
35.6
|
1.8
|
79.6
|
0.0
|
0
|
0
|
|
HAYMAN
|
2002
|
129417
|
25.4
|
0.9
|
57.4
|
0.6
|
0
|
0
|
|
GREEN CREEK
|
2002
|
4342
|
21.6
|
16.4
|
0.0
|
55.5
|
0
|
0
|
|
SPRINGER
|
2012
|
1666
|
14.7
|
0.0
|
83.9
|
0.2
|
0
|
0
|
|
RED FEATHER
|
2017
|
1582
|
13.0
|
0.4
|
79.5
|
0.0
|
0
|
0
|
|
GRACE CREEK
|
1988
|
1455
|
11.7
|
32.1
|
0.5
|
3.4
|
0
|
0
|
|
SPRING CREEK
|
2002
|
11668
|
6.5
|
0.7
|
0.0
|
45.5
|
0
|
0
|
|
HIGH MEADOWS
|
2000
|
9607
|
0.0
|
0.0
|
75.0
|
0.0
|
0
|
0
|
|
SCHOONOVER
|
2002
|
2821
|
0.0
|
0.0
|
73.2
|
0.0
|
0
|
0
|
# count number of fires with each dominant forest type
kableExtra::kable(table(fire_df$major_forest_type),
col.names = c("Dominant forest type", "# of fire events")) %>%
kableExtra::kable_styling(bootstrap_options = "striped",
full_width = F,
position = "left")
|
Dominant forest type
|
# of fire events
|
|
lodgepole
|
4
|
|
ponderosa
|
22
|
|
spruceFir
|
11
|
Explore the fire data
Map
library(leaflet)
library(rgdal)
library(geojsonio)
# read geojson file for mapping with leaflet
fire_gjson <- rgdal::readOGR(fire_filename)
## OGR data source with driver: GeoJSON
## Source: "C:\Users\tyler\Documents\Field Planning\CAREER\data\fires_stats_20210414.geojson", layer: "fires_stats_20210414"
## with 37 features
## It has 14 fields
# add a column with fire perimeter centroid point locations for map labels
fire_df$centroid_lon <- NA
fire_df$centroid_lat <- NA
fire_df$centroid_label <- NA
for (i in 1:nrow(fire_df)){
# create popup label using fire name and year
content <- paste(sep = "<br/>",
fire_df$Fire_Name[i], fire_df$Year[i])
#print(content)
fire_df$centroid_label[i] <- content
# get lon, lat coordinates of geometry centroid
fire_df$centroid_lon[i] <- st_centroid(fire_df$geometry[i])[[1]][1]
fire_df$centroid_lat[i] <- st_centroid(fire_df$geometry[i])[[1]][2]
}
# color palette
pal <- colorNumeric(
palette = "viridis",
domain = fire_gjson$Year
)
# map the fire polygons
map <- leaflet(fire_gjson) %>%
addTiles() %>%
addPolygons(stroke = FALSE, smoothFactor = 0.2, fillOpacity = 0.7,
color = ~pal(Year), group = "MTBS fires"
) %>%
addLegend("bottomright", pal = pal, values = ~Year,
title = "Fire Year",
# get rid of the comma in years (1,990 becomes 1990) for legend entries
labFormat = labelFormat(big.mark = ''),
opacity = 1
) %>%
# set base map type
addProviderTiles(providers$CartoDB.Positron, group = "Carto") %>%
addProviderTiles(providers$Stamen.Terrain, group = "Stamen")
# Icons
icons <- makeAwesomeIcon(
icon = "fire",
iconColor = "gray",
library = "fa",
markerColor = "gray"
)
# Add Fire Name popup labels to each polygon centroid
map <- map %>%
#addMarkers(lng = fire_df$centroid_lon,
addAwesomeMarkers(lng = fire_df$centroid_lon,
lat = fire_df$centroid_lat,
icon=icons,
popup = fire_df$centroid_label) %>%
addLayersControl(baseGroups = c("Carto", "Stamen"),
overlayGroups = c("MTBS fires"))
# display map
map
Fire timeline
# calculate minimum and maximum fire years
year_min <- min(fire_df$Year)
year_max <- max(fire_df$Year)
# set forest type colors for the figures
# lodgepole ponderosa spruceFir
#forest_type_colors <- c("#005a32", "#74c476", "#e5f5e0")
forest_type_colors <- c("#04663B", "#41ab5d", "#c7e9c0")
# reorder the rows based on year
fire_df <- fire_df %>%
dplyr::arrange(desc(Year))
# calculate point size - scaled relatively by fire size
tmp <- log(fire_df$Acres / max(fire_df$Acres))
tmp <- tmp + abs(min(tmp)) + 1.25
# create timeline figure
ggplot(fire_df) +
geom_point(aes(x = reorder(Fire_Name, Year), y = Year,
# color points based on dominant forest type
fill = major_forest_type),
# add a black outline around each point
color = "black", shape=21,
# size of each point, thickness of outline
size = tmp,
stroke = 0.3) +
scale_y_continuous(breaks = seq(year_min, year_max, by = 4)) +
labs(title = "Fire Event Timeline", y = "Year", x = "Fire Name\n",
fill = "Dominant forest type") +
# Put fire name on Y axis, year on X axis
coord_flip() +
# set the point fill colors using hex codes
scale_fill_manual(values = forest_type_colors,
labels = c("Lodgepole pine", "Ponderosa pine", "Spruce/Fir")) +
# set the point size in legend
guides(fill = guide_legend(override.aes = list(size=4))) +
# clean theme for plot with white background
theme_bw()

# Save figure to file
ggsave(filename = (here::here(file.path("figures",paste0(out_label,"-fire_timeline.png")))),
width = 8, height = 6)
Fire year histogram
hist_year_title <- paste("Histogram: Fire Years,", year_min, "-", year_max)
# create histogram, fire count per year.
# color the bars based on the dominant forest type of each fire event.
ggplot(fire_df, aes(x=Year, fill = major_forest_type)) +
geom_histogram(binwidth=1, color="black",
# set the outline thickness and bar transparency.
size = 0.2, alpha=0.9) +
labs(title = hist_year_title, y = "Count", fill = "Dominant forest type") +
# x axis label years from min to max year in increments of 4
scale_x_continuous(breaks = seq(year_min, year_max, by = 4)) +
# set color scale of the dominant forest types
scale_fill_manual(values= forest_type_colors,
labels = c("Lodgepole pine", "Ponderosa pine", "Spruce/Fir")) +
theme_bw()

# Save figure to file
ggsave(filename = (here::here(file.path("figures",paste0(out_label,"-fire_hist_years.png")))),
width = 7, height = 5)
Fire size histogram
# color the bars based on the dominant forest type of each fire event.
# convert from Acres to Hectares? * 0.404686
ggplot(fire_df, aes(x = Acres * 0.404686)) +
geom_histogram(aes(fill = major_forest_type), color="black",
# set the outline thickness and bar transparency.
size = 0.2, alpha=0.9) +
labs(title = "Histogram: Fire Size", y = "Count", x = "Size [Hectares]",
fill = "Dominant forest type") +
# set color scale of the dominant forest types
scale_fill_manual(values = forest_type_colors,
labels = c("Lodgepole pine", "Ponderosa pine", "Spruce/Fir")) +
theme_bw()

# Save figure to file
ggsave(filename = (here::here(file.path("figures",paste0(out_label,"-fire_hist_size.png")))),
width = 7, height = 5)
Fire size as a function of year
ggplot(fire_df, aes(x=Year, y=Acres * 0.404686 ,
log="y", label = Fire_Name)) +
geom_point(aes(# color points based on dominant forest type
fill = major_forest_type),
# add a black outline around each point
color = "black", shape=21,
# size of each point, thickness of outline
size = 4, stroke = 0.5) +
labs(title = "Fire size vs year",
y = "Hectares [log transformed]",
fill = "Dominant forest type") +
# log transform the y axis to space out the points more
scale_y_continuous(trans='log2') +
scale_x_continuous(breaks = seq(year_min, year_max, by = 4)) +
scale_fill_manual(values = forest_type_colors,
labels = c("Lodgepole pine", "Ponderosa pine", "Spruce/Fir")) +
# add Fire Name labels
geom_label_repel(aes(label = Fire_Name),
box.padding = 0.35,
point.padding = 0.5,
segment.color = 'grey50') +
theme_bw() +
# font size
theme(axis.text.y = element_text(size=16),
axis.title.y = element_text(size=16),
axis.text.x = element_text(size=16),
axis.title.x = element_text(size=16),
# move the legend to the top of the figure
legend.position="top")

# Save figure to file
ggsave(filename = (here::here(file.path("figures",paste0(out_label,"-fire_size_vs_year_scatterplot.png")))),
width = 7, height = 9)